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Summary 

Theory  behind  a  Fortran  program  (not  listed)  which  facilitates  computations 
of  transient  response  is  described.  Discrete  values  of  a  time -harmonic 
spectrum  are  used,  via  quadratic  interpolation  and  weighting  by  the 
transform  of  a  selected  time  function,  to  provide  samples  for  transient 
response  calculations  done  by  numerical  evaluation  of  a  Fourier  integral 
using  the  fast  Fourier  transform  (FFT)  algorithm.  Examples  of  motion, 
radiation  and  scattering  transient  responses  of  simple  dynamical  systems 
demonstrate  that  the  program  is  a  valuable  research  tool.  This  report  is 
an  update  of  previous  work. 
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CALCULATION  OF  TRANSIENT  RESPONSE 
FROM  TIME-HARMONIC  SPECTRUM 


By  J  H  James 
1 •  INTRODUCTION 

Transient  response  calculations  frequently  reveal  physical  features  which 
are  difficult  to  extract  from  the  corresponding  time-harmonic  spectra,  see 
for  example  James  (1).  In  order  to  facilitate  transient  sound  radiation 
and  scattering  calculations  a  Fortran  program  (2)  has  been  written  to 
produce  transient  time  histories  from  user  supplied  time-harmonic  spectra 
together  with  user  selected  time  excitations.  Numerical  examples 
demonstrated  that  the  program  is  a  valuable  tool  for  transient  far-field 
radiation  and  scattering  computations.  It  was  noted  that  while  the  program 
had  been  designed  for  calculations  of  transient  far-field  pressure  the 
amendments  required  for  it  to  be  usable  for  other  transient  response 
problems  were  trivial. 

These  amendments  have  now  been  made  and  it  is  the  purpose  of  this  report  to 
provide  an  update  to  the  theoretical  work  and  numerical  examples.  Much  of 
the  text  also  appears  in  the  previous  version  (2)  which  has  a  limited 
availability. 

2.  TRANSIENT  RESPONSE 

a.  General  Response 

The  solution  of  time-harmonic  response  problems  at  time  t  and  field-point  x 
can  be  expressed  in  the  form 

p(x,t)  =  S(x,  -u)exp(-iot)  (2.1) 

where  p(x,t)  is  the  time-harmonic  response,  which  is  the  response  due  to 
exp(-i.jt)  time  variation  of  the  excitation;  =  2~f  is  the  radian  frequency 
of  oscillation;  S(x, o)  is  the  complex  response  spectrum.  A  diagram 
illustrating  typical  responses  is  shown  in  Figure  1 . 

The  transient  response  for  separable  time  variation  Q(t)  of  the  excitation 
is  simply 

p(x,t)  =  Q('j) S(x, ij)exp( -i'jt )djj  (2.2) 

where  Q(w)  is  the  transform  of  Q(t),  viz 

Q(:->)  =  (1/2 ")/*  Q(t)exp(i_t)dt  (2.3) 

in  which  the  excitation  Q(t)  is  assumed  to  vanish  outside  the  time  interval 
t  =  0,T. 

Equation  (2.2)  can  be  evaluated  numerically  from  user  supplied  values  of 
S(x,j)  and  Q(«).  From  a  practical  point  of  view  it  is  better  written  as 

P(x>t)  =/  m  Q('o)S(x,  ,))exp(-i  u(t-tg) )  die  (2.4) 
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where  tg  is  the  transient  response  delay  factor,  introduced  in  order  to 

check  for  non-casual  behaviour  in  the  plotted  time  history;  see  Section  5 
for  numerical  examples. 


The  solution  of  time-harmonic  radiation/scattering  problems  in  the 
far-field  is  usually  expressed  in  the  form 

p(R,3,0,t)  =  A(0,tf,u))[exp(ikR)/R]exp(-iu)t)  (2.5) 

in  which  p(R,0,<p,t)  is  the  time-harmonic  pressure;  (R,6,p)  are  spherical 
co-ordinates;  k  =  ^/c  is  the  acoustic  wavenumber  where  c  is  the  sound 
velocity  in  the  fluid;  A(b,o,~)  is  the  angular  distribution  of  the 
radiation  field.  This  equation  is  often  evaluated  at  the  reference 
distance  of  one  metre,  but  many  computer  codes  simply  evaluate  A(S,o,^). 

The  transient  pressure  field  for  separable  time  variation  Q(t)  of  the 
excitation  is  simply 

p(R,0,p,t)  =  (1/R)  Q  Q(w)A(6,0,uj)exp[ik(R-ct)]da.  (2.6) 

which  is  the  same  as  equation  (2.4)  provided  S(x,^j)  =  A(9,o,ui)/R,  and 
t0  =  R/c. 

For  far-field  transient  radiation  the  standard  reference  distance  of  one 
metre  is  chosen  for  the  divisor,  R.  Note  all  responses  are  at  relative 
times  because  the  far-field  is  really  at  infinity.  This  may  cause 
conceptual  problems  because  in  certain  circumstances  causality  appears  to 
be  violated,  which  is  perhaps  not  the  case.  For  instance,  there  may  be  the 
expectation  that  the  phase  term  (R-ct)  will  cause  the  transient  pressure  to 
be  zero  before  time  t  =  R/c,  but  this  may  not  be  the  case  as  supersonic 
elastic  waves  travelling  on  the  surface  of  the  radiator  could  cause 
radiation  to  leave  a  remote  point  and  arrive  at  the  observation  point 
before  the  direct  travel  time,  see  for  example  James  (1).  Because  the 
transient  far-field  pressure  is  invariant  with  respect  to  R-ct,  it  is 
convenient  to  use  this  term  simply  as  a  device  for  shifting  the  excitation 


in  time.  Thus,  set  tr 


Rp/c  where  is  large  enough,  greater  than  a 


typical  body  dimension  say,  to  include  responses  originating  before  t  =  0. 
This  scheme  is  an  especially  valuable  facility  because  the  integration 
algorithm,  described  in  Section  4,  produces  responses  at  equidistant  times, 
starting  at  relative  time  t  =  0,  which  means  that  responses  originating 
before  t  =  0  are  lost,  unless  a  sufficiently  large  value  of  R^  is  chosen. 

c.  Far-Field  Scattering 


For  a  time  varying  plane  wave,  incident  at  the  angles  (0.,<^),  the 
insonification  is  of  the  form 

Pi(x>y.z,t)  =  /"  P^(ui)exp(-ika£)exp(-iut)duj 
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where 


(2.8) 


a  ^  =  x.  sin(e  ^)cos(j  ^)+y .  sin(H  ^)sin(.(>  ^)+z .  cos(e  . ) 

and 

P^x.y.z.t)  =  Piexp(-ikui)exp(-iwt)  (2.9) 

is  the  incident  pressure  excitation  used  in  the  time-harmonic  problem. 

Let  the  scatterer  be  surrounded  by  a  sphere  of  radius  L.  The  point  on  this 
sphere  first  met  by  the  incident  plane  wave  is  x^  =  L.  sin(0  ^)cos(p  ^) , 

y  =  L. sin(Q . ) sin($ . )  and  z  =  L.cos(o.).  At  this  point  the  incident 
s  1  s  1 

pressure  wave  is  obtained  from  equation  (2.7)  as 

Pi^Xs,ys  ’  zs’ t)  =  C,  P^w  )@xp(-ikL)exp(-iwt)dw  (2.10) 

Let  this  incident  wave  be  the  real  valued  pressure  pulse  of  finite 
duration, 

Pi(xs,ys,2s,t)  =  for  0<t<T 

=  0,  for  t<0  and  t>T 

to  give,  on  taking  the  inverse  Fourier  transform, 

P^w)  =  [exp(ikL)/27T ]/q  Q(t)exp(iut)dt 

=  exp(ikL)Q(lU) 

The  transient  pressure  scattered  to  the  far-field  is  now  synthesized  from 
the  time-harmonic  pressure  equation  (2.5)  as  the  integral 


(2.11) 

(2.12) 


p(R,e,i,,t)  =  (l/R)/“  Q(J)A(0,4,,UJ)exp[ik(L+R-ct)]4c  (2.13) 


which  is  the  same  as  equation  (2.4)  provided 
tg  =  (L+R)/c.  As  before  set  tg  =  R^/c  where 


S(x,w) 
R  is 
P 


=  A(6  ,<r  ,jj  )/R  and 
chosen  large  enough, 


>  3L  say,  to  include  responses  originating  before  t  =  0. 


3.  TRANSFORMS  OF  EXCITATIONS  IN  TIME 


a.  Impulse 

The  impulse  applied  at  t  =  0  is  defined  as 

Q(t)  =  Pi(t)  (3.1) 


where  P  is  the  magnitude  of  the  impulse.  Its  transform  is  obtained  from 
equation  (2.3)  as 

Q(io)  =  P/2  fi  (3.2) 
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b. 


The  Heaviside  excitation  applied  at  t  =  0  is  defined  as 

Q(t)  =  F0,  for  0<t<»  (3.3) 

Q(t)  =  0,  for  t<0 

Its  principal  value  transform  is 

Q(u)  =  FQi/(2™)  (3'4) 

c .  Sine  Wave 

The  sine  wave  which  is  switched  on  at  t  =  0  and  switched  off  at  t  =  T  is 
defined  as 

Q(t)  =  FgSinCu.t),  for  0<t<T  (3.5) 

Q(t)  =  0,  for  t<0  and  t>T 

the  time  for  one  cycle  being  1/fg  =  2  /«q •  Its  transform  is 

Q(J  =  (F0/2-T)[exp(i..T){i^in(.0T)-.0cos(^0T)}+.,0](.02-,2)  (3.6) 

d.  Squared  Sine  Wave 

The  squared  sine  wave  which  is  switched  on  at  t  =  0  and  off  at  t  =  T  is 
defined  as 

Q(t)  =  FQsin2(u0t),  for  0<t<T  (3.7) 

Q(t)  =  0,  for  t<0  and  t>T 

the  time  for  a  half-cycle  (one  loop)  being  1  /2f q  =  ir/ujg.  Its  transform  is 
Q  ('■>>)  =  (Fn/47riw)[exp(i.cT)-1] 

(3.S) 

_(FQ/4n:)  [exp(iuT){i.cos(a)0T)+2!o0sin(2w0T)  >-)Lo ]/ ( A-o Q2-  j2) 

e .  Attenuated  Sine  Wave 

The  attenuated  sine  wave  which  is  switched  on  at  t  =  0  and  off  at  t  =  T  is 
defined  as 

Q(t)  =  Fgexp(-at)sin(.>gt) ,  for  0<t<T  (3.9) 

Q(t)  =  0,  for  t<0  and  t>T 
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the  time  for  one  cycle  being  1/fg  =  2ir/u>g.  Its  transform  i3 
Q(“)  =  (F0/2Tr)[exp(iXT){iXsin(u0T)-to0cos(iD0T)}+u)0]/(io^-X2) 

in  which  X  =  «j+ia,  with  | a | >0 .  (3.10) 

f.  Triangular  Pulse 

The  triangular  pulse  which  is  switched  on  at  t  =  0  and  off  at  t  =  T  is 
defined  as 

Q(t)  =  (2FQ/T)t,  for  t  =  0<t<T/2 

Q(t)  =  -(2FQ/T)t+2F0>  for  T/2<t<T  (3.11) 

Q(t)  =  0,  for  t<0  and  t>T 
Its  transform  is 

QM  =  (-2FQ/2ir(D2T)  [1-exp(i,J/2)  ]2  (3.12) 

g.  Rectangular  Pulse 

The  square  pulse  which  is  switched  on  at  t  =  0  and  off  at  t  =  T  is  defined 
as 

Q(t)  =  FQ,  for  0<t<T  (3.13) 

Q(t)  =  0,  f  ;  t<0  and  t>T 
Its  transform  is 

Q(u)  =  (F0/2itiio)  [exp(iaT)-1  ]  (3.14) 

This  pulse  may  be  of  more  practical  value  than  the  impulse  and  the 
Heaviside  excitations  in  a.  and  b.  above. 

h.  Gaussian  AM  Sine  Wave 

The  Gaussian  amplitude  modulated  sine  wave  which  is  switched  on  at  t  =  0 
and  off  at  t  =  T  is  defined  as 

Q(t)  =  Fgexp(-at2)sin(u)gt) ,  for  0<t<T 

(3.15) 

Q(t)  =  0,  for  t<0  and  t>T 

Its  transform  has  no  simple  closed  form  solution,  thus  it  is  obtained  by 
numerical  integration  of  equation  (2.3). 
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The  linear  swept,  frequency  modulated,  sine  wave  which  is  switched  on  at 
t  =  0  and  off  at  t  =  T  is  defined  as 


Q(t)  =  Fgsin(at+bt2) ,  for  0<t<T 
Q(t)  =  0,  for  t<0  and  t>T 


(3.16) 


in  which  a+2bt,  the  differential  of  the  argument,  is  defined  as  the 
instantaneous  frequency.  Let  and  be  the  lower  and  upper  frequencies 

of  the  sweep,  then  a  =  and  b  =  C^-Uq)/2T.  The  transform  of  this 

excitation  has  no  simple  closed  form  solution  so  it  is  also  obtained  by 
numerical  integration. 

j-  User .Defined  Excitation 


The  Fourier  transform  of  a  user  selected  transient  excitation,  for  example 
the  amplitude  and  frequency  modulated  sine  wave  a(t)sin(b( t) ) ,  can  be 
obtained  by  numerical  quadrature  of  equation  (2.3).  This  user  pulse  is 
trivially  inserted  into  the  computer  program  (mentioned  in  Section  5)  by 
following  the  coding  of  excitations  h.  and  i. ,  provided  the  pulse 
discontinuities,  if  any,  are  at  the  limits  (0,T)  of  the  integral.  The 
coding  evaluates  integral  transforms  by  Simpson's  rule  applied  to  M  equally 
spaced  values  of  the  excitation.  When  the  excitation  has  interior 
discontinuities  the  integration  range  should  be  split  into  separate 
intervals,  bounded  by  the  discontinuities. 

4.  F.FT  EVALUATION  OF  FOURIER  INTEGRAL 

The  discretized  version  of  equation  (2.4),  which  is  of  the  form 

F(t)  =  H(Jj)exp(-iuut)d^  (4.1) 


with  H(u)=Q(io)S(x,ui)exp(ijtQ)  has  been  given  by  Cooley  et  al.  13)  as 

F(tj)  =  Gkexp(-2^ijk/N),  t.=jct,  j=0,  N-1  (4.2) 

with 

t»k  =Lp__oo  H(kou+Np5uj)  (4.3) 

N  is  an  integer  power  of  two  in  order  to  be  compatible  with  the  fast 
Fourier  transform  algorithm  (4),  lor  which  the  sampling  intervals  in 
frequency  and  in  time,  <5u  and  5t,  are  inextricably  linked  by  the  equation 

=  2ir/N«t  (4.4) 

Thus,  the  infinite  range  of  integration  has  been  divided  into  intervals  of 
length  Ni5u  to  allow  construction  of  a  quadrature  formula  which  requires  a 
single  application  of  the  FFT  algorithm  to  evaluate  the  integral  at 


N  distinct  values  of  the  time  t.  The  time  series  is  periodic,  viz  F(tj)  * 
F(tj+jj),  but  provided  the  true  time  function  decays  sufficiently  rapidly  in 

the  interval  t  =  0,N6t  it  may  adequately  be  approximated  by  discrete 
samples  derived  from  the  FFT  quadrature  formula. 


Some  trial  and  error  is  usually  necessary  before  the  constants  required  for 
numerical  integration  are  finally  selected.  The  following  suggestions  may 
be  found  useful: 


a.  First,  select  the  sampling  interval  in  frequency,  and  the 
upper  limit  of  integration,  either  by  a  guess  or  from  prior 

knowledge  of  the  spectrum  H(“).  Choose  N  (power  of  two)  and  5t  in 
accordance  with  equation  (4. A). 

b.  Stop  the  summation,  in  equation  (4.3),  at  the  integer  value 
(minimum  1)  p^  which  is  greater  or  equal  to  a ^/NSw.  Note  that  this 

may  re-define  the  upper  limit  of  integration  as  Jm"PmNlSuj. 


c.  Use  of  the  relation  H(-id)=H*(u) ,  reflecting  a  real  valued 
function  F(t),  will  halve  computation  times. 

d.  Use  not  more  than  the  first  N/2,  say,  values  of  the  time  function 
for  plotting,  as  aliasing  might  have  affected  the  values  at  the  end  of 
the  range. 

e.  If  the  time  history  plot  is  obviously  in  error,  being 
significantly  non-casual  for  instance,  then  perhaps  5w  is  too  large  or 

is  too  small. 

f.  When  the  spectrum  H(^)  is  computationally  expensive  it  can  be 
extended  with  zeros,  in  order  to  enhance  the  resolution  in  time.  Also 
the  number  of  costly  computations  of  H(")  can  be  kept  to  a  minimum  if 
an  interpolation  scheme  is  used. 

5.  NUMERICAL  EXAMPLES 

a.  General 

A  FORTRAN  program  has  been  written  to  evaluate  the  Fourier  integral 

P(x,t)  =  Q(uj)S(x,^)exp(-iuj( t-tg) )du 

numerically  from  a  user  supplied  array  of  the  discretized  spectrum  S(x,*), 
at  positive  values  of  and  a  user  selected  excitation  Q(t).  This  user 
supplied  array  need  not  be  at  equal  frequency  intervals  as  quadratic 
interpolation  is  used  to  find  values  required  for  tumerical  integration  by 
FFT.  Interpolation  is  especially  important  when  time-harmonic  responses 
are  costly.  For  instance,  if  a  spectrum  obtained  by  the  finite  element 
method  has  sharp  peaks  it  may  be  desirable  to  form  a  composite  spectrum 
from  a  crude  resolution  over  the  entire  frequency  range  combined  with 
several  highly  resolved  regions  where  peaks  are  dominant.  The  program 
automatically  pads  the  user  supplied  spectrum  with  zeros  if  its  highest 
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frequency  is  less  than  the  maximum  frequency  ljJm=PraN'5u  in  the  FFT  quadrature 

formula.  This  useful  feature  allows  selection  of  a  time  resolution  which 
is  small  enough  to  obtain  smooth  plots  of  the  transient  response,  without 
the  added  expense  of  unnecessary  high  frequency  computations  of  S(w). 

For  the  numerical  examples  the  following  table  summarizes  the  sampling 
constants  used,  frequency  values  being  in  hertz: 


I  Spectrum  Constants 

1 _  _  _ 

FFT  Constants 

j  Fig;  NF 

FMIN 

FMAX 

f 

N 

t 

f 

f 

m 

'  2a  ! 1001 

0.25 

1000.0 

1.0 

i 

8192 

0.000500 

0.244 

2000 

2b  201 

0.25 

1000.0 

5.0 

i 

8192 

0.000500 

0.244 

8000 

13  1001 

0.25 

1000.0 

1.0 

i 

16384 

0.000125 

0.244 

8000 

1  4  11001 

0.25 

1000.0 

1  .0 

i 

16384 

0.000125 

0.244 

2000 

>5  !  501 

1  .00 

1000.0 

2.0 

i 

2048 

0.000250 

1.953 

4000 

!6  1001 

10.0 

1000.0 

.99 

i 

8192 

0.000200 

0.610 

_ 

5000 

where  for  the  time-harmonic  spectrum  constants,  NF  is  the  number  of 
frequency  values  ranging  from  FMIN  to  FMAX  with  a  resolution  of  5f;  and  for 
the  FFT  constants  f  is  the  upper  frequency  limit  of  integration.  Padding 

is  brought  into  effect  in  all  examples  because  the  upper  frequency  limit  in 
the  FFT  constants  is  greater  than  the  upper  frequency  limit  in  the  spectrum 
constants.  At  the  upper  frequency  limit  of  each  spectrum  the  product 
Q(w)S(w)  is  small  enough  to  allow  truncation  of  the  Fourier  integral 
without  any  significant  loss  of  accuracy.  Computational  times  for  all 
examples  are  trivial,  when  compared  to  the  times  needed  to  evaluate  the 
spectra  S(u) . 

When  the  material  is  steel  the  following  SI  units  are  used;  Young's  modulus 
19.5  x  1 0 10 ,  Poisson's  ratio  0.29,  shear  modulus  7.558  x  1010  and  density 
7700:  water  has  density  1000  and  sound  velocity  1500, 

b-  Impulse  Response  of  Mass-Soring 

In  the  first  example  the  grounded,  heavily  damped,  mass  and  spring  system 
shown  in  Figure  2  is  subject  to  a  unit  impulse  force.  The  constants  in 
SI  units  are  mass  1000,  stiffness  1.0  x  10e  and  hysteretic  loss  factor  0.1. 
The  displacement  response  spectrum  S(u)  was  calculated  by  a  program,  called 
COUPLE  (5),  which  evaluates  the  response  of  a  complex  dynamical  system  by 
combining  dynamic  stiffness  matrices  of  simple  dynamical  systems.  The 
spectrum  S(w)  is  smooth  except  in  the  region  of  the  mass/spring  resonant 
frequency  at  50.3  Hz.  The  spectrum  Q(w)  of  the  excitation  is  flat. 

The  impulse  response  is  shown  in  Figure  2a,  the  time  delay  being  selected 
to  start  the  response  at  t  =  80  ms.  It  is  in  excellent  agreement  with 
Milne  (6)  who  shows  mathematically  and  demonstrates  numerically  that  the 
very  small  non-causal  response  just  before  t  =  80  ms  is  due  to  the 
selection  of  hysteretic  damping.  A  plot  (not  shown  here)  has  also  been 
obtained  in  which  the  damping  is  viscous  rather  than  hysteretic:  the 
non-causal  response,  as  expected,  is  not  present. 
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Figure  2b  shows  non-causal  response  behaviour  introduced  by  under-sampling 
the  response  spectrum  S(u).  In  this  case  the  frequency  spacing  in  S(uj)  is 
too  large  for  accurate  interpolation.  The  response  is  acceptable  for 
practical  computations,  but  it  could  be  improved  dramatically  by  inserting 
about  four  additional  S(u)  values  centred  on  the  50.3  Hz  natural  frequency. 

C-  Asia!  response  of.  FM  Excited  Beam 

In  the  second  example  the  uniform  beam  shown  in  Figure  3  is  excited,  in  the 
axial  direction,  by  a  linear  FM  sine  wave  force  of  unit  amplitude  whose 
frequency  sweeps  from  50  to  100  Hz  in  20  ms;  see  excitation  in  Figure  1. 

The  steel  beam  is  of  length  20,  area  cross-section  0.012996  and  has  an 
hysteretic  loss-factor  of  0.01.  The  acceleration  response  spectrum  S(u)  at 
each  of  the  beam  tips  was  calculated  by  COUPLE.  They  vary  rapidly,  with 
prominent  resonances  at  126  Hz  spacing.  The  spectrum  Q(u)  of  the 
excitation  drops  rapidly  above  400  Hz. 

The  axial  acceleration  responses  at  the  ends  of  the  beam  are  shown  in 
Figure  3,  the  time  delay  being  selected  to  start  the  responses  at 
t  =  10  ms.  The  responses  are  substantially  causal  before  t  =  10  ms.  The 
response  of  the  left-hand-side  tip  in  Figure  3a  resembles  the  excitation 
initially,  and  then  extended  ringing  is  evident,  the  beating  being  due  to 
interaction  among  the  first  three  resonances.  Compare  with  the  response  of 
the  right-hand-side  tip  in  Figure  3b:  the  plots  are  almost  identical  after 
relative  time  40  ms,  when  allowance  is  made  for  a  delay  of  4  ms  for  the 
wave  to  reach  the  far  end  of  the  beam  at  the  longitudinal  wave  speed  of 
5032  ms  l;  levels  in  the  former  plot  being  just  a  fraction  higher  because 
of  the  effect  of  damping. 

d-  Impact  Response  of  Isolation  System 

In  the  third  example  the  simple  isolation  system  shown  in  Figure  4  is 
excited  at  its  base  by  a  unit  triangular  pulse  force  which  is  on  for  10  ms. 
The  machine  is  modelled  as  a  steel  Timoshenko  beam  in  bending  with  length 
1.0,  area  cross-section  0.012996,  moment  of  area  in  bending  1.407  x  10'5, 
hysteretic  loss-factor  0.01  and  shear  correction  factor  0.822.  The 
isolators  are  modelled  as  simple  springs,  which  do  not  transmit  bending 
moments,  of  stiffness  0.333  x  108  and  loss-factor  0.01.  The  foundation  is 
modelled  as  an  infinite  steel  plate  of  thickness  0.01  and  loss-factor  0.02. 
The  vertical  acceleration  responses  S(u)  were  calculated  by  COUPLE.  They 
vary  slowly  except  for  a  prominent  peak  at  564  Hz  which  is  due  to  the  first 
bending  resonance  of  the  beam.  There  is  no  sharp  mass-spring  resonance 
because  the  machine  is  mounted  on  a  low  impedance,  but  highly  resistive, 
foundation.  The  spectrum  Q(u>)  of  the  excitation  drops  smoothly  to  zero  at 
200  Hz  and  its  first  side  lobe  at  300  Hz  is  27  dB  down. 

The  vertical  acceleration  responses  of  the  tips  of  the  beam  are  shown  in 
Figure  4,  the  time  delay  being  selected  to  start  the  responses  at 
t  =  80  ms.  The  response  of  the  left-hand-side  tip  in  Figure  4a  resembles 
the  excitation,  with  small  perturbations,  for  the  first  10  ms  and  then 
extended  ringing  is  evident.  Compare  with  the  response  of  the  right-hand- 
side  tip  in  Figure  4b:  the  initial  response  is  of  opposite  phase  while  the 
extended  ringing  is  of  the  same  phase.  These  observations  are  consistent 
with  the  early  response  being  dominated  by  the  heave  and  pitch  modes  of  the 
beam  and  the  later  response  being  due  to  ringing  at  the  fundamental  beam 
resonant  frequency.  Interaction  through  the  water  loaded  plate  is  probably 
small. 
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e. 


Radiation  and  Scattering  of  Shell 


The  fourth  example  is  of  far-field  radiation  from  and  scattering  by  the 
thin  elastic  steel  shell  shown  in  Figure  5.  The  excitation  is  a  transient 
force  or  plane  wave,  of  unit  amplitude,  which  has  the  form  of  a  300  Hz  sine 
wave  which  is  on  for  one  cycle.  Shell  constants  are  radius  1.0,  thickness 
0.01  and  hysteretic  loss-factor  0.01.  Radiation  and  scattering  spectra 
S(w)  were  calculated  from  closed-form  theoretical  expressions  (7,8).  They 
vary  rapidly  in  the  mid-frequency  regime  due  to  resonance  effects.  The 
spectrum  Q(u)  of  the  excitation  has  a  broad  lobe  near  to  250  Hz  and  its 
first  side  lobe  near  725  Hz  is  18  dB  down. 

Figure  5  shows  plots  of  the  far-field  pressure  transients.  Small  ripples 
at  times  before  15  ms  are  caused  by  premature  truncation  of  the  spectra 
S(ui).  For  the  radiation  transient  the  initial  pressure  has  the  same  phase 
as  the  excitation,  probably  due  to  the  force  being  applied  directly  to  the 
water.  For  the  scattering  transient  the  initial  pressure  has  the  opposite 
phase,  probably  due  to  the  incident  wave  initially  seeing  a  pressure 
release  shell.  Extended  ringing  caused  by  resonances  is  very  marked  when 
it  is  recalled  that  the  excitation  lasts  for  only  3.3  ms. 

These  plots  also  appear  in  the  earlier  publication  (2)  and  there  they  are 
accompanied  by  plots  of  spectra  S(cu)  derived  from  both  closed-form  and 
finite  element  solutions. 

f .  Radiation  from  Water-Filled  Pine 

The  final  example  is  of  far-field  airborne  sound  radiation  from  the 
water-filled  steel  pipe  shown  in  Figure  6.  The  excitation  is  a  transient 
point  source,  of  unit  amplitude,  which  has  the  form  of  a  250  Hz  sine  wave 
which  is  switched  on  for  either  one  cycle  or  forty  cycles.  The  pipe  has 
radius  0.10,  thickness  0.005  and  hysteretic  loss-factor  0.01.  The 
radiation  spectrum  S(w)  was  calculated  from  closed-form  theoretical 
expressions  (9).  It  varies  smoothly  except  for  sharp  resonances  at  225  and 
685  Hz.  The  spectrum  Q(ui)  of  the  excitation  has  its  main  lobe  near  250  Hz. 
Figure  6  shows  plots  of  the  far-fieid  pressure  transients.  Small  ripples 
at  times  before  the  main  response  starts  are  caused  by  premature  truncation 
of  the  spectra  S(oj).  In  Figure  6a  ths  excitation  is  a  single  cycle  sine 
wave;  the  initial  response  closely  resembles  the  excitation,  and  the  later 
response  is  dominated  by  slowly  decaying  ringing  of  the  225  Hz  resonance. 

In  Figure  6b  the  rapid  transition  to  a  near  steady  state  is  shown.  The 
beating,  slowly  decreasing  with  time  around  the  time-harmonic  level  of 
±  0.44,  is  caused  by  interaction  between  the  continuous  excitation  and  the 
decaying  resonance  transient.  When  the  excitation  is  switched  off  after 
40  cycles  the  amplitude  drops  rapidly,  leaving  decaying  ringing  due  to  the 
resonance  transient. 

6.  CONCLUDING  REMARKS 

Transient  response  (motion,  radiation,  scattering,  etc)  of  a  dynamical 
system  can  be  expressed  as  a  Fourier  transform  over  frequency  of  the 
time-harmonic  pressure  weighted  by  the  transform  of  the  time  varying  part 
of  the  excitation,  assumed  separable  in  space  and  time.  This  transform  is 
evaluated  numerically  by  an  application  of  the  FFT  algorithm.  The  main 
features  that  have  arisen  from  the  work  reported  herein  are  as  follows: 
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a.  A  FORTRAN  program  which  evaluates  a  transient  response  from  a 
user  supplied  array  of  time-harmonic  responses,  which  need  not  be  at 
equidistant  frequency  intervals  as  interpolation  is  used  to  obtain  the 
values  required  for  numerical  integration.  The  array  of  pressures  may 
be  extended  automatically  with  zeros  in  order  to  improve  resolution  in 
time.  Computer  time  is  trivial. 

b.  A  number  of  basic  transient  excitations  which  can  be  selected  by 
the  user  of  the  program.  It  is  not  difficult  to  insert  additional 
excitations  into  the  computer  program,  even  if  closed-form  expressions 
are  not  available  for  their  transforms,  as  a  scheme  for  numerical 
integration  is  built  into  the  program. 

c.  Numerical  examples  of  transient  responses  which  provide 
guidelines  for  the  program' s  use  and  also  demonstrate  that  it  is  a 
powerful  research  tool. 

Computations  of  transient  response  of  linear  systems  are  usually  done 
either  by  'time  marching'  on  the  original  time-dependent  differential 
equations  (or  matrices)  or  by  Fourier  transform  of  the  time-harmonic 
solutions,  as  herein.  For  general  response  calculations  the  latter  method 
would  appear  to  have  an  advantage  not  only  because  computer  programs  for 
time-harmonic  problems  are  widely  available  but  also  because  the  appearance 
of  non-causal  responses  are  an  indication  of  errors,  usually  due  to 
under-sampling  or  premature  truncation  of  the  time-harmonic  spectrum  S(w). 
Time  marching  methods  implicitly  assume  c?usality. 

Future  extensions  to  the  program  should  include  a  limited  capability  for 
filtering  the  transient  response,  as  is  usually  the  case  in  physical 
measurement.  A  causal  filter,  with  frequency  response  spectrum  F(w),  is 
easily  inserted  in  the  computer  program  as  the  transient  response  is  now 
the  Fourier  transform  of  F(w)Q(u))S(w)  ■ 
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AMPLITUDE  IN  MICROMETRES 


small  nan-causal  dip  due  to 
assumed  hysteretic  damping 
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FIG.  2  Impulse  response  of  damped  mass-spring  system: 

(a)  overdetermined  spectrum  sampling  gives  excellent 
results  and  (b)  underdetermined  spectrum  sampling  gives 
a  non-causal  response. 


ACCELERATION  IN  MILLIMETRES/SECOND  ACCELERATION  IN  MILLIMETRES/SECONt? 


FIG.  4  Vertical  acceleration  response  of  isolation 
system  excited  by  triangular  pulse  on  for  10ms: 

(a)  at  near  beam  tip  and  (b)  at  far  beam  tip. 


PRESSURE  IN  MILLIPASCALS  PRESSURE  IN  MILLIP/.SCAL9 


RELATIVE  TIME  IN  MILLISECONDS 


FIG.  6  Far-field  airborne  pressure  in  mi22ipasca2s 
ref.  lm  of  water-filled  pipe  excited  by  250Hz  sine 
wave  point  source  on  for  (a)  1  cyc2e  (b)  40  cycles. 
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